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Abstract 

This work summarizes the most important developments in the con- 
struction and apphcation of the Dirac-Moshinsky oscillator (DMO) with 
which the author has come in contact. The literature on the subject is 
voluminous, mostly because of the avenues that exact solvability opens 
towards our understanding of relativistic quantum mechanics. Here we 
make an effort to present the subject in chronological order and also in 
increasing degree of complexity of its parts. We start our discussion with 
the seminal paper by Moshinsky and Szczepaniak and the immediate im- 
plications stemming from it. Then we analyze the extensions of this model 
to many particles. The one-particle DMO is revisited in the light of the 
Jaynes-Cummings model in quantum optics and exactly solvable exten- 
sions are presented. Applications and implementations in hexagonal lat- 
tices are given, with a particular emphasis in the emulation of graphene 
in electromagnetic billiards. 
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1 Introduction 

The harmonic oscillator is the paradigm of integrability and solvability with 
applications to many branches of physics. As it was written by Moshinsky [T] 
"...A complete analysis of the subject would require an encyclopedia, within 
which one of the volumes could be (our) book" . Is it possible to promote all 
these features to a relativistic quantum-mechanical model? This is the question 
that Moshinsky and Szczepaniak answered in a seminal paper more than twenty 
years ago. Today we can see how this idea has been exploited in several ways 
by using the Dirac-Moshinsky oscillator (DMO) as a way to understand better 
the mathematical structure of solvable Dirac equations. But beyond the math- 
ematical developments surrounding this system, in these lecture notes we would 
like to emphasize that some applications can be found in areas of physics of 
current interest, such as the study of electrons in two dimensional materials (for 
example, graphene) and the interaction of atoms with electromagnetic fields in 
cavities (the Jaynes-Cummings model). The notes are divided in four sections. 
In the first section we give a detailed introduction to the subject, covering sym- 
metries, Lorentz covariance and algebraic solvability. In section two, we review 
the many body theory for the Dirac equation (in first quantization) and the key 
points of the spectral structure of these systems. We continue with the formula- 
tion of the Dirac oscillator as an interaction between Dirac particles and a brief 
mention to hadronic spectroscopy is made. In section three we present solvable 
extensions to the single particle DMO in the context of isospin fields and con- 
tinue with the formulation of an exact mapping of such extensions to quantum 
optical cavities. Finally, in section four we deal with tight binding lattices, two 
dimensional systems and the effective Dirac equations appearing in materials 
such as graphene and Boron Nitride. In the same section we develop the same 
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idea in the context of electromagnetic billiards and a deformation method is 
proposed, leading to a realization of a DM0 in one and two dimensions. 

2 The Dirac-Moshinsky oscillator for one parti- 
cle 

2.1 The Dirac Oscillator as proposed by Moshinsky and 
Szczepaniak 

Our purpose is to review the construction of an interaction for relativistic sys- 
tems (particles) producing bound states for arbitrarily high energies with ana- 
lytically solvable spectrum. Lorentz invariance is crucial. This was achieved by 
Moshinsky and Szczepaniak (1989) [2] with further generalizations to describe 
interacting particles [S] through Poincare invariant equations. 

A naive approach to the problem is to propose a one-particle relativistic 
equation in the form 

{c^h^A + ni^c^ + Knoj'^r^)(f) = (1) 
with the trivial result that the energies become 

E"^ -m^c^ = 2Loh{n+^) (2) 

However, the Lorentz invariance of the problem is not clear in this simple picture. 
It is also necessary to find a first order equation in time (as Dirac originally 
proposed through his equation [6|) for a good application to initial condition 
problems, as it is the case for hamiltonian systems in quantum mechanics. 

Interactions which are linear in the coordinate were introduced in but 
Moshinsky and Szczepaniak introduced and solved a Dirac equation with a 
hamiltonian of the form 

= ca • (p ± iujm(3r) + m(?f3 (3) 

where p = — i?iV and the Dirac matrices are given by /3 = 7'^, a' = /37*, 
i = 1, 2, 3. The 7's, in turn, are given in the usual representation 
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The double sign in the frequency oj written in indicates that similar results 
can be obtained independently of this choice. In this framework, both coordi- 
nate and momentum operators must appear in linear form in order to preserve 
integrability: a clear indication of phase space symmetry. The symmetry Lie 
algebra of this system was investigated in |22[ . and the corresponding genera- 
tors are now represented in the algebra of Dirac matrices. The corresponding 
group decomposes naturally into 0(4) (compact component representing a non- 
relativistic oscillator) and 0(3, 1) (non-compact component representing states 
with infinite degeneracy). To see how these two types of degeneracies appear, let 
us analyze the stationary solutions of the Dirac equation with the hamiltonian 
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2.2 Stationary solutions 

The stationary form of our equatfon H"^ = E"^ has bispinor sohitions of the 
form 



? ) (5) 



satisfying 
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+ TOw^r^ + mc^ + 3uj + 4-L • S W-a = -B>2 (7) 

TO ft 



where we have used the customary definition S — ^Ticr. The wavefunctions are 
given in terms of the isotropic harmonic oscillator states with total number of 
quanta N = 2n + I and orbital angular momentum I. Such states are coupled 
to the spin ^ as we now indicate: 



^i=ANji\N{l,^)jm) (8) 

2c 

and Aj^ji is a normalization constant. The energies result in 



ip2 = — (E + mc^) • (p — iTOwr)?/'! (9) 
n 



and we write the wavefunctions associated to the positive and negative energies 
in the form 

The completeness of these eigenfunctions pT|) has been proved in [7] as a 
straightforward exercise. See the figure [1] for an explanation of the two possibil- 
ities of the spectrum according to the parity of the orbital angular momentum 
I. Here it is worth to mention that these solutions constitute a way to write 
a propagator in spectral form, and that the wavefunctions themselves can be 
computed through the exact expression of the Dirac oscillator Green's function, 
obtained in [51 [T7]. 



2.2.1 Non-relativistic limit 

Using our previous relations, it is easy to see that 



{E^ - m^c'^)iji = 
(c2(p2 + w^m^r^) - Shcumc^ - ^J^rnc^'i- ' s) V'l (12) 

for which the relativistic energy given by e = i? — rac^ <^ mc^ leads to 
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Squares of energy levels as a function of n,v,j 

E 



(5.1) 


(4.2) 


(3.3) 






(6.0) 


(7,1) 


(8,2) 


(9.3) 


(10.4) 


(3.1) 


(2.2) 








(4.0) 


(5.1) 


(6,2) 


(7.3) 


(6.4) 


(1.1) 










(2.0) 


(3.1) 


(4.2) 


(5.3) 


(6.4) 


(0.0) 


(1,1) 


(2.2) 


(3.3) 


(4.4) 



v=2 



n=2 



Figure 1: Structure of the spectrum for the 3D Dirac-Moshinsky oscillator as 
proposed by Quesne and Moshinsky [H]. The eigenvalue e (see (fT5)) ) is shown 
as a function of the total angular momentum and as a function of parity in 
alternating rows. The quantum number n = \{N — j + ^) parameterizes the 
infinitely degenerate states. The number ly = ^{N + j — 3/2) gives the states of 
finite degeneracy. The nomenclature (a, b) corresponds to a = 2N + l,b = I for 
infinite degeneracy and a = 2v + I ,h ~ I for finite degeneracy 



ei>i = {hho - ~ 2jL ■ yj^ (13) 

where Hho is the usual harmonic oscillator hamiltonian. This shows that the 
non-relativistic limit reduces to the oscillator without its rest energy (the con- 
stant term —^hoj in (jl3p ) and with a strong spin-orbit coupling. The infinite 
degeneracy does not disappear, but the negative energy solutions decouple from 
small components of the spinors as expected. This leaves us with ipi as our 
non-relativistic states. 

2.3 The Dirac-Moshinsky oscillator in Lorentz covariant 
form 

In relativistic problems, it is always important to formulate everything in Lorentz 
covariant form. This ensures that the solutions obtained in a particular frame 
of reference (such as the expressions obtained before) are valid in other inertial 
frames under the appropriate Lorentz transformations. For our purposes and in 
what follows, it is convenient to rescale all quantities such that our units give 
h = c = 1. Furthermore, let us work with uj — I and leave the mass m as 
the only free scale of our system. The Lorentz covariant wave equation for the 
DM0 can be given as 
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(7" [Pm - *^±M^-7''] + m) * = (14) 

where 7'^ are the Dirac matrices as defined in ^ and the perpendicular projec- 
tion of vectors is given by 

r±f_c = r^- {r''u^)u^, (15) 

the vector Ui, being a time-hke four vector such that (uj/) = (1,0,0,0) for some 
inertial frame. In such a frame of reference, (jl4p can be written as 

H^ = ^- (16) 

with H given by It is tempting to regard the 4- vector potential —ir±^Ui,'-f'' 
as a minimal coupling with a gauge field; however, we must warn the reader 
that the matrix l3 not only precludes this possibility, but also ensures that such 
a minimal substitution is not a "pure gauge" interaction, therefore giving non- 
trivial results. We can give a physical meaning to the vector u^, as we show in 
the next section in connection with anomalous coupling. 

2.3.1 Pauli coupling 

The anomalous (Pauli) coupling is a way to introduce interactions with an 
external field (say, a magnetic field B) producing terms of the form S • B in the 
hamiltonian. It is not our purpose to delve into the nature of such a coupling, but 
we may emphasize that it provides a possibility of preserving gauge invariance 
other than the usual minimal coupling. Using the definition of the spin tensor 
S^i^ — (l/4){7^, 7iy}, the Dirac equation for our DM0 can be written as 



[7^p^' + m + Sf^^F^^ = (17) 

with the choice F^'^ — u^r'^ — u'^r^. The meaning of the external field F can 
be found by noting that 



d^F^- = -u^ (18) 

i.e. the vector u'^ can be interpreted as a current. The Maxwell equations for a 
field given by the tensor F suggest that a constant current given by the r.h.s. of 
(|18p would produce a Dirac oscillator by means of anomalous couplings. Finally, 
in the frame of reference (1,0,0,0) our current gives a uniform charge density 
filling the space. 

2.3.2 The supersymmetric formulation and its extensions 

A supersymmetric algebra |20) which has the squared hamiltonian as its center 
can be identified as the responsible for the infinite degeneracy of the DM0. We 
will revisit this point in further sections using a different notation. For now, 
let us recall what has been done in the context of supersymmetry. When a 
non-abelian vector potential is used to produce a Dirac oscillator (i.e. p i— > 
p -I- i/3A(r)), one can prove the relations 
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{Qa,Qb}^Sab{H^ -1), [Qa,H^]=0 (19) 

with 

0>=(."a *=(.;.a -'V')- ™ 

This structure reveals that more than one choice for the vector potential A(r) 
allows analytical solvability: For a general expression of the form a = p + 
iG{r)r, the radial function G{r) may lead to a harmonic oscillator or a Coulomb 
problem, both of them with additional centrifugal barriers. This is related to 
the factorization method devised by Infeld and Hull [9^, as it was noted in [50] 
in connection with the radial equation resulting from the substitution of a in 
the Dirac equation. One has the radial equation 

and the choices G{r) = a/r+b, a' /r+b'r are possible, leaving the supersymmetry 
algebra intact. In order to break the degeneracies one may try several tricks. 
In particular, the introduction of interactions depending explicitly of the total 
angular momentum is a way of breaking such degeneracies by hand. It is also 
evident that this approach is attached to the dimensionality of the problem since 
the radial equation has been used to propose the corresponding extensions. 
In the following, we shall use an alternative approach to understand infinite 
degeneracies in connection with dimensionality (2 or 3 dimensions). It will 
result that the one-dimensional Dirac oscillator admits a superalgebra similar 
to the one given above, but its degeneracy (if any) is strictly finite. 

2.4 Hilbert space and algebraic structure 

Here we introduce a notation and some concepts which cast the DM0 as a bi- 
linear form in bosonic and fermionic operators. This will prove useful in the 
discussion of invariants and spectral properties. The Lorentz group is locally 
isomorphic to SU{2) x SU*{2). The Hilbert space of our problem is therefore 
L2{C) X S3X S3, where 6*3 is the space of normalized complex vectors of two en- 
tries (Pauli spinors) . Obviously, each 5*3 is a three-sphere [TU] . In the following, 
our hamiltonian will be given by 

H = a{p + il3r) + mP (22) 
with the following representation of the Dirac matrices 




For reasons that will become apparent in further sections, we may refer to this 
representation as quantum-optical. With this notation we may introduce the 
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concept of *— spin through the vector S^, whose z-projection eigenvalues account 
for big and small components of spinors. Upon rotations, this projection also 
gives solutions with positive and negative energies. 

^+^(0 o')=^+®^2, S_ = (S+)t, E3 = /3 (24) 
The Hamiltonian can be written in algebraic form as 

= I]+S-a+I]_S-a1'+mE3, (25) 

The dependence of H on ladder operators makes evident the fact that the fol- 
lowing operators are invariant: / = a''' • a+ 5S3, I' = (a- tr)'''(a- a) + With 
the integrals of the motion given by a combination of fermionic and bosonic 
operators, we can obtain the sohitions of the eigenvahie problem as follows. 

Two states with angular momentum j and satisfying the eigenvalue equation 
I\ ) = (2n + j — 1)1 ) are given by 

|</)i) = |n, (i - 1/2, l/2)i,m,-)|-), \cf>2) = |n - 1, {j + 1/2, l/2)j, m^)\+). (26) 

Another pair of states with the same angular momentum j but with eigenvalue 
7| ) = {2n + 3)\ )is 

1.^3) = |n, (i + 1/2, l/2)i,m,-)|-), |<^4) = |n - 1, {j - 1/2, l/2)j, m,)|+). (27) 

The 2x2 blocks of H obtained from these states can be obtained easily. Here 
we give such blocks 



H(,,2n + ,) = { ^5=1!^ \/25r+7)y (29) 
V ^/2{n + j) m J 

Prom the eigenvalue equation applied to these subspaces, we obtain the well 
known energies E'^ = ni? + 2(n + j) and E"^ = + 2ri, which correspond 
to the expressions we have found before. Infinite and finite degeneracies come 
from these two blocks respectively. Let us now go further and write similar 
expressions for low-dimensional Dirac oscillators. 

2.4.1 Boson- Fermion algebra for 1 + 1 and 2 + 1 dimensions 

The discussion on the algebraic structure above can be implemented directly 
in 1 + 1 and 2 + 1 space-times. For this we have to find the boson (harmonic 
oscillator) and fermion (*-spin) operators which parallel our previous discussion. 
For the 1 + 1 case we define a, in terms of the position x and the momentum 
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p in the standard form. For the 2 + 1 case, it is useful to define the following 
chiral creation and annihilation operators (subindex r for right and I for left) 



Ur = ai + ia2, ai = ai — ia2 = (a^)* (30) 
with the properties 

[cir, ai] = [ttr, (ai)*] ^ 0, [ur, 4] = [ai,aj] = 4. (31) 
The low dimensional hamiltonians are 

=ai{p + iPx) + mP, (32) 

with ai = — CTi, P = as and 

F(2) = J2 aiiPi + iPn) + mp, (33) 

i=l,2 

with the low dimensional Dirac matrices chosen as ai = —(72, a2 = —ai,P = a^. 
These hamiltonians can be written in algebraic form as 

fl'(^) =a+a + <T_at + mcT3 (34) 

H^^^ = a+ttr + cr-a\. + mas (35) 

Both of them have a 2 x 2 structure: The spin is abscnit in one spatial dimension 
and a± corresponds to *— spin, while in two dimensions as also generates the 
U (1) spin leading to the total angular momentum L3 + ^as- The solvability can 
be viewed again as a consequence of the existence of invariants. In this case we 
have 

/(i)=ata+ia3 (36) 



= ar-al + i(T3, Js = Urol - Ojaj + ^as (37) 

The two dimensional case exhibits some peculiarities. The conservation of angu- 
lar momentum J3 comes from the combination of cr and in iJ '^^^ , together with 
the absence of ai,aj. This absence is also responsible for the infinite degeneracy 
of all levels. On the other hand, the three dimensional example is manifestly 
invariant under rotations due to its dependence on S • a and S • a'^ and its infinite 
degeneracy comes from the infinitely degenerate operator (cr • a)((T • a)^'. 

Let us summarize the material of this section. We have learned that the 
intcgrability of the harmonic oscillator can be implemented in the context of 
the Dirac equation by recognizing that coordinates and momenta should lie on 
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an equal footing - the essence of phase space symmetry. The eigenfunctions 
and energies were given exphcitly. A Lorentz covariant equation with a Dirac 
oscillator potential could be written and interpreted in terms of anomalous 
coupling and a constant external current. The infinitely degenerate part of the 
spectrum could be understood either in terms of a supersymmetric algebra (in 
the 3 + 1 dimensional case) or as a consequence of the non-compact part of the 
symmetry Lie algebra (unitary representations are infinite dimensional). We 
went further and gave a description of the Dirac oscillator in terms of fermionic 
and bosonic ladder operators (the operators a, a'' and <7±), showing thus the 
existence of integrals of the motion in a more transparent way. It was also shown 
that the degeneracies in the 3 + 1, 2 + 1 and 1 + 1 dimensional examples obey a 
different pattern; in three and two dimensions the parity plays an important role 
(absence of j and absence of left chiral operators, respectively), while the one- 
dimensional DM0 cannot have infinite degeneracy in despite of the existence of 
a supersymmetric algebra (one spatial degree of freedom is insufficient). 

3 The many body Dirac equation 

The success of Moshinsky's work related to the harmonic oscillator of arbitrary 
particles and dimensions is due to the fact that the results provided a good basis 
to solve variational problems in bound composite systems [1]. This was imple- 
mented in composite models describing atomic nuclei. The idea is to extend 
this success to relativistic quantum mechanics, with the obvious application to 
many-particle systems where high energies are involved. Many of these examples 
can be found in the context of hadron physics, where "relativized" models have 
been proposed [291 [30]. However, we need a model which allows the integrabil- 
ity and solvability we are seeking for, in order to understand the structure of 
multiparticle relativistic formulations, rather than just fitting the results to ex- 
perimental data. To this end, we start here with the many body Dirac equation 
as proposed in |5] , followed by the study given by Moshinsky [TH] , [H] regarding 
the positive part of the spectrum and degeneracies in a general framework (the 
Foldy-Wouthuysen transformation [13], [IS], [Ij). Then, we present the two 
and three particle Dirac oscillators as proposed again in the list of works [5] . 

3.1 Poincare invariance of the many body problem 

We review the generalization of the Dirac equation for a system of many particles 
(carefully treated in pi). The main idea is to mimic the treatment of many 
particles in non-relativistic quantum mechanics as the direct product of operator 
spaces for particle l,2,..,n. The main equation is defined such that, in the frame 
of reference where the center of mass is at rest, we recover a hamiltonian of the 
form 



where Hi is the Dirac hamiltonian of the i-th particle. The potential V is 
assumed to be independent of the center of mass. Such an equation is 



N 




(38) 
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N 



ip^O. (39) 



The relative coordinates and the time-hke relative coordinates are given respec- 
tively by 



— Xfii ^±fj, ^ "^^ij (40) 

The meaning of the time-like vector defining our preferred frame of reference 
is obvious, as the hamiltonian stands for the energy at the center of mass with 
four vector P^. We must use the time-like unit vector in the form 

= (-P,P-)-i/2p^. (41) 
For convenience we have defined the matrices 



N 



r = nT>M> r, = (7>^)-ir. (42) 

r=l 

Taking and iJ = P° in one recovers (pS)) . 



3.2 The Foldy-Wouthuysen transformation 

The problem of positive and negative energies in the Dirac equation for arbitrary 
potentials appeared from the very beginning and it was treated systematically 
by Foldy and Wouthuysen [T3], [H]. However, in most of cases such a treatment 
can be carried out only approximately. See, for example, |;15^ for a detailed 
review of the subject. Remarkably, the Dirac oscillator is one of the examples 
(together with the free case) in which the corresponding transformation can be 
carried out analytically (211 . There exists a unitary operator which transforms 
the Dirac hamiltonian into a diagonal operator in spinorial components. In our 
algebraic language, the transformation finds the basis in which the z component 
of *-spin gives the positive and negative energies of the system. The idea is to 
express the hamiltonian in terms of its even part (diagonal matrices) and odd 
part (anti-diagonal matrices) and find a hermitian operator S such that 



Hfw = e'^Hoe^'^ = even. (43) 
For the free particle one has iS = P{a ■ p)0, tan(20a • p) = a • p and 



Hfw = Np^ + m2 (44) 

For the three dimensional DM0 we use the definition a • tt as the kinetic energy 
of the DM0. One has now the relations iS — /3{a ■ ti)9, tan(20a • tt) — a-Ti and 



Hfw = /3^p2 + r2 + (3 + 2L • cr)/3 + to2_ (45) 
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In the following we give a more detailed treatment dealing with arbitrary poten- 
tials, first for one particle and then for many particles. This will be useful in our 
interpretation of a many particle Dirac equation based on the direct product of 
particle spaces. 

3.3 The many body Foldy-Wouthuysen transformation 

With the aim of characterizing the spectrum of a multibody system with inter- 
actions, we seek for an expansion of H in terms of inverse powers of the rest 
mass. Such an expansion should allow the identification of positive and negative 
energies of the model. For one particle in a potential V, we have 



H = + £ + V, O =cxp, £ =m/3. (46) 
We apply a unitary operator U = exp(iS') exp(iS") exp(iS"'), 

2m 2m 2m 

0' = A[a.p,y], 0" = zi|-P>!, H' = UHU^ (47) 

Expanding up to l/{mass)^ in the kinetic energy, l/{mass)'^ in the potential, 
we have 



H' = H + V, H = /3lm + ^-^ 
(p X E) - (E X p) 



(48) 



with E = —W, S = -j^ot X a. For two particles we define the corresponding 
matrices as 



ai=a(8)l, q;2 = 1 (8) a (49) 
/?i=/3(8l, /32 = l®/3. (50) 

The hamiltonian is H = Hi + H2 + y(ri,r2). Applying successively Ui = 
exp(zS'i) and 1/2 = exp{iS2) one gets 

C/2C/i/f(C/2C/i)^ = ^i+ ^2 + + higher order (51) 
In the general case with n particles, one has 

JV 

H = J2Hi + V{ri,...,rN) (52) 



JV 
i=l 
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with 



t= l,2,---n 

At the end, we have an expression which shows the first relativistic corrections 
to the kinetic energy, the spin-orbit couplings and the Darwin term [23j . But the 
most important feature of the result is that the positive energies can be extracted 
immediately by selecting the positive eigenvalues of all the Pt multiplying the 
kinetic energy. 



3.3.1 The cockroach nest: Extraordinary infinite degeneracy 

Before completing our task of generalizing the Dirac oscillator to many parti- 
cles, it is important to understand first the types of degeneracies involved in a 
multiparticle Dirac equation. Here we show a very simple example. For com- 
muting Dirac hamiltonians one expects that the total FW transformation can 
be decomposed into individual factors corresponding to each hamiltonian. Ac- 
cording to our definition of the multiparticle FW transformation, the free case 
gives 



Hfw =^Pi\lpi + (55) 

1=1 

where it becomes evident that the energies are now added with 'wrong' signs due 
to the (3 matrices. This means that the transformation to even hamiltonians 
contains both particle and anti-particle solutions without a correction of the 
signs in front of their kinetic energies. Specifically, for two particles of equal mass 
described by an observer at the center of mass, only the relative momentum p 
appears. One of the corresponding energy eigenvalues has the form ^/p^ + m? — 
\f(^-p^^^^rrr? = for any p. Moreover, when the relative momentum vanishes 
the rest energy of the system becomes instead of the usual value of 2m. This 
result seems to be unphysical. Therefore, one has to project the final result onto 
the purely positive component, otherwise one would obtain an extraordinary 
infinite degeneracy [24j. 



3.3.2 Application to the two body problem 

As a point of comparison and before dealing with integrable problems, let us 
consider now a system of two particles with an interaction given by a quadratic 
potential. The system is not integrable. The transformed hamiltonian can be 
approximated by 



1 

4m^ 



Si + S2 



(p X E) - (E X p) 



(56) 



with p the magnitude of the relative momentum. The potential is so far arbi- 
trary. We may propose a quadratic interaction V = ^ma;^(ri — r2)^. 
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In the center of mass frame, the choice of positive energy components reduces 
the hamihonian to 



H'=(2m + 3f) + (P- + ^ + fs.L)-J^ (57) 
\ am J \ni 4 4to / 4m'' 

with r, p the relative coordinate and momentum, respectively. The spectrum of 
the problem is found by diagonalizing the matrix with elements given by 



{n'l, (^iij5;j,m|i/'|n/, (^iij^,j,m) 

= {2m+ ^+i^(2n + I + l) + f^bXj + 1) - Kl + 1) - sis + 1)] ]6„ 

1 /,n';'i^4| 



{n'l'\p^\nl) 

where we use two-particle harmonic oscillator states with spin, i.e. 



(58) 



|nZ, Q ^; J, m >^ ^ < ^M, Sa\jm) \nlfi) | Q Sa) (59) 

We take N < Nmax to get a finite matrix. 

As an application, one can describe the mass spectrum of binary systems 
such as bottomonium or charmonium. It is very important to comment on the 
flavor of the wavefunctions: According to the theory of particles composed by 
quarks (Hadrons), one has to use appropriate wavefunctions containing the in- 
formation of quark flavor and color, the interaction being flavor-blind |27| . In 
other words, our potential is permutationally invariant. Moreover, for approxi- 
mately degenerate quark masses (for example u and d quarks) the states must lie 
in isospin doublets and should be properly (anti)symmetrized according to the 
representations of the permutation group [3S]. However, for quarkonia we only 
have one pair made of quark-antiquark, and the process is trivial. It is there- 
fore sufficient to consider products of wavefunctions of the form Flavor xSpinor 
which are symmetric, since the colorless feature of the composite demands an 
antisymmetric color part. In the following, all symmetric and antisymmetric 
solutions of the eigenvalue problem related to our hamiltonians are considered. 

It is possible to introduce quartic corrections to the potential in order to 

4 4 

obtain more realistic spectra V' — — '""^ '' . The FW transformation of such 
a term yields next order corrections, therefore we neglect them. The coupling 
constants and the rest mass are taken as adjustable parameters. They are fitted 
to experimental data [31] using least dispersion. See the figures. 

3.3.3 Application to the three body problem 

The hamiltonian is now 



'-E/3^(™* + ^-£!3)+i^«Hp.xE.-E.xp.H^V?T/+y (60) 
with a flavor-blind potential 
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C^o") Bottomonium System 

m=4.5 Gev, (j=0.34 Gev, a=0.92 



T(3S) 



T(2S) 



T(1S) 



Figure 2: Energy comparison. Solid line: experiment. Dashed line: theory. In 
, J denotes the total angular momentum (also referred to as total spin) and 
P the parity of the corresponding energy state 



V 



Mw2 



(ri - + (r2 - ra)^ + (ra - ri)^ 



(61) 



Now we use the definition of Jacobi coordinates in order to separate the contri- 
bution from the center of mass of the system and the relative coordinates. The 
harmonic oscillator for n particles with hamiltonian 



2 " 



(62) 



can be decoupled into n — 1 oscillators by using the Jacobi coordinates in the 
form 



{Ps)j = [s(s + l)]-i/2^((pt),-(p«+i),),s = l,...,n-l, 

t=i 

n 

{Pn)j = n-y'^{pt)j. (63) 



Considering that the center of mass is at rest, we obtain 



1 .2 ,1 



1 1 



H' = Mc''+(]pi + ^pi) (— + —)+- P2 + ^P12. 

V4 12 / V"^i "^2/ 3rn3 V12 \mi m2. 
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Charmonium System 

m=1.4 Gev, (j=0.51 Gev, a=0.05 



■/'(2S) 



J/*'(1S) 



X^l (IP) 



X,, (IP) 



Figure 3: Energy comparison. Solid line: experiment. Dashed line: theory. In 
, J denotes the total angular momentum (also referred to as total spin) and 
P the parity of the corresponding energy state 
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1 1 1 1 T. 
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3\/3 



3m: 



S3 • Lo 



2^3 



(64) 



where L12 = r 1 x p2 + 1 O 2 and pi2 = Pi • P2- 
The spectrum is obtained by diagonalizing 



{n[,l[,n'^,l'^,L'; ( ^^^T'^S';j'm'\H'\n,,h,n2,l2,L; ("iiVl^; jm) (65) 



where the states are 



\ni,li,n2,l2,L; [ ^^)^^'5';jm) = 
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Sigma Baryons 




Figure 4: Dashed: Theory. Sohd: ExperimentaL To achieve a better agreement 
with experiment, we have considered a shght variation of the parameters w, M 
as functions of the integrals of the motion (see figures [2j [3] for the meaning 
of this nomenclature) . The corresponding values can be found in table [U with 
the particular feature that the effective total rest mass is close to 1.2Gev 



^ < Lfi,Sa\jm)\ni,li,n2,h,Lfi)\(--jT-S(T > (66) 

The matrix elements are computed by means of Racah algebra. We take N„iax = 
3. Again, this application does not include other degrees of freedom such as 
particle flavor. However, this does not preclude the use of these results to obtain 
a part of the spectrum for a three quark system in which constitutive masses 
are considered, with the result that the members of our composite system are 
distinguishable particles (with very different masses and broken degeneracy). 
As it is evident, the hamiltonian we use in this case is not permutationally 
invariant, in despite of the fact that the potential enjoys of such a property. 

To achieve a better agreement with experimental data, we may introduce a 
mass and a frequency which depend on the integrals of the motion. We include 
a comparison with the spectra of S particles (strange baryons) . In summary, we 
have shown that the Poincare invariant equation with arbitrary inter-particle 
potentials can be treated in the quasi- relativistic approximation by means of the 
Foldy-Wouthuysen transformation, with the possibility of computing spectra 
of the transformed hamiltonian. From the computational point of view, the 
process is not necessarily simple and in despite of the few parameters that can 
be used to fit energy levels, our understanding of the system is now beyond the 
symmetry principles underlying the Dirac-Moshinsky oscillator. In the following 
we describe how to construct the DM0 hamiltonian for more than one particle 
and analyze the corresponding solutions. 
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JP uj (Mev/?i) M/1.2 Gev 



1 + 

2 

1- 
2 

3 + 
2 

3- 
2 

5 + 
2 

5- 
2 



96 

184 

187 

179 

137 

137 



1.00 
1.00 
1.27 
0.93 
1.11 
1.03 



Table 1 : Table of parameters 



3.4 The two-particle Dirac oscillator 

Now we proceed to generalize the DM0 to two particles. There is more than 
one generalization which gives a solvable two-particle problem, as can be seen in 
[T] , [S] . Here we concentrate in one possibility for the interacting potential. For 
simplicity, let us set the masses of the particles as unity. It is also convenient to 
restore our frequency lo in this section, as we want to analyze the spectrum in 
terms of the coupling. After all, one can always go back to the former units by 
replacing cj ^ Lofi/mc^. The hamiltonian and the Poincare invariant equation 
are, respectively 



H 



[oLi - da) • (p - i-rB) + Pi + h 



(67) 



^ r, (7^^^^ - iwxl^^r) + 1) 



* = 



(68) 



The interaction matrix B is chosen here as ,9i/^275i752 with 751 — 75 1, 
similarly for 752. In this case, the hamiltonian given above admits an expression 
which is quite similar to that of our previous algebraic analysis. With the 
appropriate definitions of *-spin operators S3, E§ for particles 1 and 2 

and the bosonic operators a, for the relative coordinate, one has 



- h.c. + 



(69) 



showing clearly that the structure leading to infinite degeneracy is still present 
through Si • a. The infinite degeneracy of the cockroach nest manifests itself by 
the "wrong" addition of energies. One could repeat the treatment given before 
in terms of invariants, which are still easy to identify. However, the original 
work of Moshinsky (described with detail in [1]) did not rely on this possibility 
and proceeded in the direction of decomposing the spinors and iterating the 
resulting equations connecting them (this corresponds, implicitly, to compute 
the fourth power of the hamiltonian) . The resulting spectrum can be found as a 
function of the total angular momentum, total spin and total oscillator quanta: 
E = ±EN,s,3,m with 
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-ojN,0 



(2VT 

EN,s,j,n, = { 2 ^1+c^(7V + 2 



for s = 0,P = (-)J 
for s = 1,P = 
[ 2v^l +a;(7V + 1),0 for s = 1, P = -(-p 



(70) 



The wavefunctions are known for all cases indicated before. Let a±,b±, c±,± be 
numerical coefficients; then we have the wavefunctions 



V "^22 j 



_ 1 / a+ + a_ 



\N{j,0)jm) 



(71) 



valid for s — 0. Whenever s — 1 and P — (—)-' , we have 



V'li 

'ip22 



+ 6- 
- 6_ 



|iV(j, l)jm) 



(72) 



For s = 1, P = -(-)J the result is 



-011 
i>22 



72 V - c+ 
%/2 



c 



|iV(7 + l,l)jm) 
|iV(j-l,l)jm) 



(73) 



where the coefficients c±±,a±, b± are determined by the secular equations aris- 
ing from the Schroedinger equation for the relativistic hamiltonian. Taking into 
account (j67p . the stationary equation yields the complementary components of 
the wavefunction V'2i,V'i2- These wavefunctions can be found in terms of Racah 
coefficients in the appendix of |32) . 

3.5 The three-particle Dirac oscillator 

For this case, we follow again the book by Moshinsky [T| and recognize that the 
center of mass can be eliminated from the outset by proposing the following 
Poincare invariant equation and hamiltonian: 



,r) + 1] U = 



(74) 



(75) 



where the primed observables denote the operators for a particle of index s 
after subtracting the corresponding observable for the center of mass (either the 
momentum or the coordinate operators). Here, the matrix B in the interaction 
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is chosen as /3i (g) . . . (g) /3„ , and we may choose m particular n = 3. The spectrum 
is obtained by combining the equations for some of the spinor components of 
the wavefunction and by noting that the total number of quanta of two Fock 
states is conserved (corresponding to the oscillator states of the two remaining 
Jacobi coordinates). The wavefunctions are 



*4 



/ V'lii \ 






i>122 






i>212 




■0211 


\ i>221 J 




\ ■0222 / 



and they satisfy 



with 



0, O = MDZ^M^ - 



(76) 



(77) 



diag {E -"i.E + 1,E + 1,E + 1), 



(78) 



diag [E -l,E ~l,E -l,E + ?,) 



(79) 



/ S3 • a^, 
S2 • a'2 
Si • a'l 

V 



S2-a^ 
Ss-aJ, 


Si • a'l 



51 • a'l 


S3-a(, 

52 • a'2 



\ 
Si • a'l 

S2-a^ 

S3 • a^ ; 



with the operators without the center of mass defined as 
a^ = a^ - -(ai + a2 + a3). 



(80) 



(81) 



In this treatment we recognize again the pattern of *-spin provided by the 
components cx |±), leading to a hamiltonian of the form H — S+M + 
h.c. + E3 + I]| + Unfortunately, this does not lead to interpretations which 
are similar to our previous examples. The reason one has to diagonalize the 
operator O instead of using algebraic properties of the hamiltonian is related 
to integr ability: In this case we have 8 invariant operators given by P^,3,N, 
while the total number of degrees of freedom (without taking into account the 
spin of the particles) is 9. 

The application of this problem given by Moshinsky et al. [T] [S] consisted on 
the calculation of the spectrum of masses of nucleons, together with a compar- 
ison with their experimental masses. There, the light quarks u, d were treated 
as identical particles. Using the irreducible representations of the permuta- 
tion group, suitable wavefunctions of the form Flavor x Spinor were used in 
the computation of energies and eigenf unctions. The application went as far as 
computing a form factor for the proton by using the information of the resulting 



20 



N 






n.1 


n2 


h 


h 


P 


L 


J 























+ 





S 


1 


1 











1 







1 


\1-S\<J<1 + S 


1 





1 











1 




1 


|1-5'| < J< 1 + 


2 


2 





1 











+ 
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< L < 2 


\L-S\<J<L + S 
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2 











2 
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\2-S\<J<2 + S 


2 





2 











2 


+ 


2 


\2-S\<J<2 + S 



Table 2: Table of states for Nmax = 2 



wave function of the ground state. It is not our purpose to repeat such a prowess 
here, but instead let us compute the spectrum of this system without any other 
degrees of freedom than the ones provided from the outset. Our intention is to 
analyze the scaling properties when the coupling oj is varied from small to large 
values. Using the states 

|ni,/i,n2,Z2(i);~(r)i(5);JM) = 

[ihlnih) X {f2\n2l2)]L X 

one can find the matrix elements of O . The resulting matrices are finite for 
each number of total quanta. We restrict to N = 0,1,2. The wavefunctions 
can be finally obtained by finding the null vectors of the matrix {O ) for each 
energy. The complementary components are obtained as before, i.e. by using 
the original stationary equation. See the table of states. 




-4-2 2 

E 



Figure 5: Energies for lu = 0.03. The eigenvalues are distributed in four groups 
around the values —3, —1, 1, 3. The states coming from the cockroach nest lie 
around 1,-1 

In summary, the two and three particle Dirac oscillators can be solved. In the 
two-particle case, the features of the spectrum could be identified straightfor- 
wardly, given the simplicity of its hamiltonian. Moreover, we could show that 
such a hamiltonian obeys the algebraic scheme proposed for the one-particle 



(l|^)x(2|i) 



X (3I-) 



(82) 



TA/T 
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Figure 6: Spectrum for w = 0.1 and N = 2. As the interaction increases, the 
levels become more spaced. 



E 

Figure 7: Spectrum for w = 1 and N = 2. Here we see a typical spectrum for 
which the interaction does not allow a quasi relativistic expansion 

case, exhaustively analyzed in section 1. The three-particle DM0 constitutes a 
more challenging example, since the number of integrals of the motion does not 
match the total number of degrees of freedom. However, the numerical diago- 
nalization can be done without much effort. The resulting spectra showed levels 
grouped around energies 1,-1, which can be regarded as a consequence of the 
cockroach nest. In the following, we comment on the n particle case, where we 
expect similar results. 

3.6 One dimensional n particles 

Let us discuss briefly the infinite degeneracy present in this model. For our 

purposes, we may eliminate the spin of the particles by restricting ourselves to 
one-dimensional space. The only degrees of freedom remaining in our simplifi- 
cation are given in terms of the annihilation operators without center of mass 
and the *-spin operators associated to each particle. We argue that the kinetic 
part of the DM0 hamiltonian 
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n 

H ^{1 + B)'^ (j\a[ + h.c. + mass 



(83) 



is infinitely degenerate. To see this, we may apply eigenstates \s\) of (j\ to 
(i/ — mass)^. Take states of the form -0 — |s})...|s^) x |A^i)...|A^„_i) where each 
\Ni) is a Fock state with the only condition that ^ Ni = constant. Clearly, the 
definition of our operators a[ gives Y^a'^ = Q and any of the states proposed 
above for which all the s} are equal gives a vanishing kinetic term, regardless of 
the total number of quanta. The mass term thus not touch the oscillator parts, 
therefore the resulting matrix elements using these states lead to an infinitely 
degenerate spectrum. The cockroach nest makes itself present for an arbitrary 
number of interacting particles. Its elimination is not a trivial task, in despite of 
our careful choice of observables. As a final question and in view of the results 
presented in this subsection, one may ask whether a system of one-dimensional 
Dirac particles can parallel the Calogero model |33], which is one of the most 
general integrable models in the non-relativistic realm. 

To end this section let us quote Moshinsky regarding the applications of 
his work on relativistic oscillators: "We conclude by stressing that we have 
made a calculation using a harmonic oscillator picture with a single parameter 
(frequency) and it is as good or as bad as many more complicated ones that 
start from QCD or that use many more parameters." 

4 Exactly solvable extensions 

In this section we present extensions of the one-particle DM0 which allow solv- 
ability and connect our relativistic systems with implementations in Quantum 
Optics, via the Jaynes-Cummings model [31] for one and two atoms. The solv- 
able extensions that we propose can be motivated entirely in the framework of 
non-local relativistic potentials [33], [33], [37], [35], [35]. The key point in the 
introduction of such potentials is the presence of isospin: An internal degree 
of freedom of our fermion which we shall couple to an external field [19] , 02] . 
Interestingly, many of the problems related to infinite degeneracies of the DM0 
(discouraging the use of its wavefunctions as a basis for more complicated prob- 
lems) can be removed in a natural way, preserving the simplicity of the model. 
We suggest the reader to follow references HT]. 

Consider a hermitian operator of the form $(r,p) as the potential to be 
introduced in the total hamiltonian. One has H^'^^ = h'^^ -\- $, with h'^'' given 
by the dimensional Dirac-Moshinsky oscillator treated previously. On phys- 
ical grounds, this corresponds to a bound fermion perturbed by a momentum- 
dependent potential. We introduce also an internal group for this field, for 
example the SU{2) associated to isospin or as the gauge group of a non-abelian 
field. Let us denote the corresponding Pauh operators for isospin by Ti,T2,T3, 
with the usual definitions for the ladder operators T± . The simplest expression 
that we can use is a linear one in a and has the form 

$ = (T+S -a + T^S -a^ +7T3) (84) 
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where we now use 7 to denote a coupling constant. In fact, one may consider any 
potential of the form $ = F(r+S-a+T_S-at+7T3) where F is a function which 
admits a power expansion. Clearly, [A^+ iTg, $] = . A suitable group of states 
can be used to evaluate the 4x4 blocks of H. We describe this procedure by 
restricting ourselves to the linear case (|84| for simplicity. The lower dimensional 
examples follow the same pattern 

(A + rraB) (r+a + T_a+ +7T3) (85) 

i?^^^ — Cr+ar + (J-al + 771(73 + 

{A + asB) {T+ar + T^a]. + 7T3) (86) 

H^^^ = S+S • a + I]_S • a^ + 771S3 + 
(A + S3B) (T+S -a + T^S -at +7T3) . (87) 

With these extensions, it is evident that the new invariants for one, two and 
three dimensions are 

/(I) = a) a + ^(73 + ^Tg (88) 
/(^) = aral + ^0-3, J3 + ir3 = arol - aia\ + ^0-3 + ^^T^ (89) 



=at .a+-S3 + -r3, J=atxa + S. (90) 
4.1 Analytical Spectrum 

Now we compute the eigenstates of H'^^^ . We evaluate the 4x4 matrix H{N, j) = 

|0f) = |77,(j + l/2,l/2)j-, 777,)h)s|-)T (91) 

\ct>^) = \n,{j-l/2,l/2)j,m,)\-)^\+)T 
10^) - 177 - 1, (j - 1/2, l/2)j-,7n,)|+)s|-)T 
K) = |77-l,(j + l/2,l/2)j,7n,)|+)s|+)T 

where 77 is the oscillator radial number, j is the total angular momentum and 
777j its projection in the z axis. These are eigenstates of /^^^ with eigenvalue 
TV = 277 + j - 1/2. 

The resulting 4x4 blocks of H with elements H{N,j)ki = ) are 
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I -m,-{ A-B)-f {A-B)^2{n + j) -^2{n + j) 
{A - B)^2{n + i) -~m+{A-B)-i 

-v/2(n + j) m-(A + S)7 (A + B)^/2^ 

V (^ + B)\/2^ m+{A + B)jJ 

and the energies can be obtained explicitly for each of these blocks using the 
formula for the roots of a quartic polynomial. The infinite degeneracy is now 
broken, since one cannot reduce H{N) to smaller blocks where only n appears. 
The exception to this occurs when A = B = 0, recovering the DM0 without 
additional external fields. 

4.2 Lorentz invariant form and Pauli coupling revisited 

With the aid of a vector we can introduce more interactions in a covariant 
way. A non-local, non-abelian field tensor J^^" = X^^^i^i-^f be intro- 
duced in the equation by means of the Pauli coupling. We propose 

^^"^^eA^-^P^Ar^p (92) 
_^P-^^p.Ap^^p^^ (93) 

^r=0, (94) 

for which the Dirac equation reads 

[7^p^ + TO + 5^,F^'^+B5^,^n^-0, (95) 

This type of fields have been introduced with the purpose of describing a 
finite characteristic length (due to non-locality) and also as a way to prevent 
divergences in perturbation theory. The nature of such fields can be elucidated 
by inserting our J^^y in the corresponding non-local field equations [3S], [55] , 
[ID]. Using 



T''" =u^irlTi+plT2)- (96) 



one has 



^^'^ = j(K, i?"] - i/) + [B^, B^]. (97) 
The gauge potential and the current can be obtained in the form 

Bp. = Ufj.i^rf^r'^Ti + r„p'^T2) Bilinear in p, r. (98) 

r = *[pp,j-n + [i3M,-^n (99) 



.2' 

-u'^Ti + + trilinear terms in p,r 



= -u'^T, +p1 + ( -W,r^<} - K,rl}r^ ) (100) 
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We do not elaborate further on these points, since our aim here is to simply note 
that our momentum-dependent potentials admit a treatment which is parallel 
to that given in section 1. We refer the reader to [55] for detailed derivations. 
The eigenvalues for the one-dimensional extension are given in the figure [SI 



Energy levels as a function of coupling 




i 



Coupling (Discrete values^ 
Levels without external f eld 



Figure 8: Spectrum of our exactly solvable extension. The vanishing coupling 
shows the eigenvalues of the Dirac oscillator. Degeneracies are lifted and level 
spacing increases. 

In this section we have shown that the introduction of extra degrees of free- 
dom motivated by particle physics (e.g. isospin) can be used to propose exactly 
solvable models. Furthermore, such extensions can be rewritten in Lorentz in- 
variant form and the coupling of the external fields to our DM0 was shown 
to parallel the Pauli coupling treated in section 1. Also, the removal of infi- 
nite degeneracy was possible by extending the space (rather than introducing 
a j dependence in the interactions) . Although we have computed the spectrum 
analytically, it is desirable to understand the dynamics of this system. Such a 
problem can be analyzed by methods dealing with entanglement of the different 
observables of our system: spin, *-spin, isospin and oscillator operators. 

4.3 Quantum Optics 

The remarkable analogy between the Dirac oscillator and the Jaynes-Cummings 
hamiltonian has been pointed out before [H] , [33] , [li] with the aim of producing 
such a system in a quantum optical experiment. The structure of our extended 
hamiltonian shows that our model can be mapped to a Jaynes-Cummings hamil- 
tonian of two atoms (of two levels each) inside an electromagnetic cavity. If the 
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dimensions of such a cavity are properly chosen, the eigenmodes of the quantized 
electromagnatic field will be sufficiently separated in frequency, with the possi- 
bility of coupling our atoms to only one boson operator. The one-dimensional 
example and the double Jaynes-Cummings model coincide: 



H = a+a + (T^a) + mcra 4- T+a + T^a) + 7T3 (101) 

where we have to identify cr, T with the operators for the two-level atoms 1 and 
2. The operator a is now the annihilation operator of the electromagnetic field 
mode. Spin-spin interactions can be introduced as well. 



4.3.1 Dynamical application: Entanglement and Decoherence 

The origin of entanglement and decoherence measures is related to quantum 
information and quantum computation |46| . However, such quantities can be 
defined and computed in such a simple way that they can be used to analyze 
the dynamical features of general systems involving several degrees of freedom. 
Here we shall take advantage of this situation and proceed to define a partition 
of the system A + B. 

We take a pure state density operator p = \^{t)){ilj{t)\ of the entire system 
and compute purity P and entropy S of the Dirac oscillator subsystem. 

P{t) = TrN,. ((Tr,p(t))') 

S{t) = -TrN,.(Tr,p(t)Log(Tr,p(t))), 

(102) 

where TrN.o- is the trace with respect to oscillator and *-spin degrees of freedom, 
while Tr^- is the trace with respect to isospin. Let us analyze the one dimensional 
case for simplicity. The integral of the motion is 



JW =ata+i^3 + ir3 (103) 



We use the eigenstates of /^^^ 



|0?) = |n + 2>|--> m) = \n+l)\-+) 

|0«) = |n + l)| + -) \^2) = \n)\++) (104) 



H = 



Ho 











Hi 











H2 



\ 



(105) 



where Hn is a 4 x 4 block. 

Now we analyze the entanglement of a Dirac oscillator with the external 
field, the initial state is chosen a.s ip — Xn ® Xj 
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Ix) = l/\/2(cos6'|+) +sin6i|-)) 
and Xn is a solution of the unperturbed Dirac oscillator 



(106) 



\Xn)^Ai+^n)\+)+Ai-^\n+l)\-) (107) 

We use the exact energies and wavefunctions to compute P{t), S{t) (purity and 
entropy). Other initial conditions can be used in the context of Quantum Optics, 
for example in cases where the initial state is prepared as a product of the two 
atoms |±)|±). In that case, the external field induces entanglement between 
such degrees of freedom, although the atoms do not interact directly but only 
through the cavity. However, this side of the analogy will not be discussed here 
gH [H]. here. 



Purity Purity 




Time^ natural units. Time, natural units. 



Figure 9: A resonant effect around 7 = m. The field produces entanglement in 
a regime where the energy is nearly the rest mass. We have used purity as the 
simplest way to characterize the entangled state. 



Entropy 



Entropy 





Time natural units 



» so 
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Figure 10: A resonant effect around 7 = to. The field produces entanglement 
in a regime where the energy is nearly the rest mass. Here we show the entropy 
by way of comparison. The structure obtained in figure [S] is also present in this 
result 
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We have learned in this section that the path of algebraic solvability leads 
to many possibilities regarding relativistic wave equations with additional de- 
grees of freedom such as isopin and their ad hoc realization in quantum optical 
experiments where atomic levels can be used to emulate certain observables. 

In our dynamical study, we have seen that the perturbation of a DM0 with 
fields of arbitrary intensity admits an analytical treatment, with the conclusion 
that the rest mass of the system responds to the stimulus via entanglement. 
Our toy model suggests that particle creation and maximal entanglement are 
related. On the other hand, this can be interpreted as a resonant effect in the 
Quantum Optics analogy. 

In what touches the experimental realization of this model, the following 
setup can be proposed with quite general parameters. We require two atoms of 
different species trapped in an electromagnetic cavity. The coupling constant 7 
defined in relation with m (fixed as unity), can be adjusted by placing the atoms 
in regions with different field intensities. Moreover, a large distance between the 
atoms is needed in order to ignore direct interaction terms. The evolution of 
entanglement studied in this section can be realized by preparing the initial 
state. In principle, one could prepare it by trapping one atom in the cavity and 
measuring the energy of the total system (this corresponds to a dressed state of 
the Jaynes-Cummings model with one atom). After this is achieved, the second 
atom can be introduced in the trap. This setup thus emulates our model and 
makes it experimentally accessible. 

5 Emulating a Dirac-Moshinsky oscillator in Elec- 
tromagnetic Billiards 

In this section we give an account of our recent findings related to hexagonal 
lattices and the emulation of Dirac equations. In recent years, there has been 
an explosion of papers (for instance, |50| and references cited therein) related 
to the experimental observation of true monolayers of carbon obtained from 
graphite: graphene. The technique, known as micromechanical cleavage, takes 
advantage of the property that graphite is composed of weakly interacting layers 
of carbon atoms and such layers can be removed and analyzed individually in 
atomic microscopes (graphene flakes). Together with the many possibilities for 
the practical applications of such materials, there is an additional feature which 
is of particular interest to our subject: Relativistic quantum mechanics. 

The band theory of graphite was studied by Wallace in the 1950's [49] . 
There it was shown that the dispersion relation of electrons propagating in a 
hexagonal lattice becomes linear at the edges of the Brillouin zone - a hint for 
a relativistic energy formula, although the slope is given by the Fermi velocity 
instead of the speed of light. At the corners of the Brillouin zone, one could 
find conical energy surfaces and, according to the quantum field theoretical 
approach proposed by Semenoff [5T] , one could also find a pseudospin from the 
decomposition of the hexagonal lattice into two triangular sublattices. Since 
it is the hexagonal structure what is essential to this analogy, one could also 
try to emulate the same behavior by propagating waves in periodic arrays of 
resonators - whose resonances should play the role of the atomic orbitals in 
a material. This happy analogy meets the existing technology of microwave 
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cavities, originally used in the context of chaotic billiards. Here we review the 
corresponding theoretical treatment in a detailed manner and go further by 
proposing a model of deformations giving rise to an effective wave equation 
given by the Dirac-Moshinsky oscillator. 



Carbon atoms 




electronic structure 



Figure 11: Deformation of a graphene sheet, induced by the proximity of 
Lithium atoms (blue spheres). A hexagonal cell is shown in red. This plot 
can be obtained using Density Functional calculations to determine the equilib- 
rium positions of the atomic centers. See, for example, Seligman and Jalbout 

SHI). 



5.1 One-dimensional Dirac equation 

The situation described above can be modelled in a simple way by a Schroedinger 
equation with a potential consisting of deep wells, each of them located at a 
lattice point. The specific shape of atomic wave functions is irrelevant, as long 
as we know how the overlaps (interactions) decay as a function of the distance 
between resonators. For practical purposes, such decay can be regarded as 
exponential, which follows from considering a lattice of constant potential wells. 
As an additional remark, such potentials should be deep enough such that only 
one level (or isolated resonance) contributes to the dynamics. 

A lattice consisting of two periodic sublattices is considered. They have the 
same period and are denoted as type A and type B. Each sublattice point can 
be labeled by an integer n according to its position on the line, i.e.Xn- The 
energy of the single level to be considered in the well is denoted by a for type 
A and /3 for type B. The state corresponding to a particle in site n of lattice 
A is denoted by \n)A and the corresponding localized wave function is given by 
^a{x — Xn) — {x\n)A- The same applies for B. The probability amplitude A (or 
overlap) between nearest neighbors is taken as a real constant. 



The hamiltonian of a tight-binding chain can be cast in terms of Pauli ma- 
trices (Ta , (T_|_ = ai + ia2 , U- =0"^ by defining 
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5.53 6.50 6.70 6.8C 6.Q0 

frequency/GHz 



Figure 12: Density plot of the coupling between two resonators. The exponential 
decay of the coupling is demonstrated by noting that two resonators constitute 
a two level system for which E ^ Eq ± A. As the distance increases, the level 
splitting A tends to zero and the two peaks (orange paths) merge into a single 
peak exponentially fast. Courtesy of U. Kuhl. 



n : 



A A 
A A 



V 



(109) 



and setting Af = (a - /3)/2, Eq = {a + (3)/2. We have 



H = Eo + asM + a+U + a^U'' (110) 

This is a general structure which explains the appearance of pseudospin. 

It is left to show that there is a region where the spectrum is linear (Dirac) . 
The spectrum is computed by squaring H. 

{H - Eo f ^ + im^ (111) 

Bloch's theorem enters in the form 

U(l>k = A(l + e'2->fc)0,., nnVfc - A^ll + e'2-Afe|2^^ (^^2) 
The energies and eigenfunctions of H are 



E{k) ^Ea±JA'^\l + e»2^^fe|2 + AP 



(113) 
(114) 
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Figure 13: Configuration of potential wells (or resonators) on a chain. We use 
yellow and black discs to recall that in principle the resonators can be of different 
type, a) The periodic case, b) General deformation, c) Dimer deformation. 



Around points where the inter-band distance is minimal, we have the usual 
relativistic formula 



E{k) = Eo± v'A2k2 + M2, (115) 

The amplitudes are proportional to the overlap between neighboring sites 
and decay exponentially as a function of the separation distance between res- 
onators, i.e. 



Ae-'*"/^, (116) 



where dn + X is the separation distance between resonators of type A and B 
in the n-th position. When d„ = 0, the periodic configuration is recovered. 
The length A has been introduced for phenomenological reasons: The decay law 
might be given by a multipole law, but we fit it to an exponential decay by 
adjusting A. 

With all this, it is natural to expect a modification in the operators 11,11^'. 
We use a, and impose [a, a^ = uj = constant (The limit a; = recovers 
Bloch's theorem). One finds the conditions 



A„,„ = A, A2+i_„+2 - A2 „+i = oj (117) 
Therefore the distance formula for the resonators is 
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Figure 14: Resonators in a one dimensional lattice. The plot above gives a 
representation of resonators as a function of the x-coordinate, while the plot 
below shows an ideahzation of the corresponding potential (wells) and the wave 
functions of resonances. These functions may leak outside the wells. 



/ A2 \ 

dn = Alog -r^ , < n < Umax (118) 

with n.nax = [l4rl]- 

Finally, we have the hamiltonian 

H ^ En + a'iM + a+a + (T^a'^ (119) 
with energies and wave functions 



E{n) ^Eq± y/ojn + AP, 0>n> A^/lo (120) 



V' - ±(E{7i)-Ea)-M ; 



5.2 Two-dimensional Dirac equation 

The concepts given in the last section are now extended to produce an emulation 
of graphene. We shall use the same algebraic strategy to derive spectra and a 
possible extension through deformations, namely the two dimensional Dirac- 
Moshinsky oscillator. 

5.2.1 The free case in 2D 

We start with the definition of the vectors which generate our hexagonal lat- 
tice (see figure [T71) . It is divided in two triangular sublattices, one of them 
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Ground State, Real Amplitude 






Site Number 



Figure 15: Ground state as a function of site number. The ground state wave- 
function is obtained by multiplying the values given in the ordinate by the 
individual resonant wavefunctions. These are considered to be highly peaked at 
each site. The signs alternate from site to site. The envelope is approximately 
gaussian (nodes are absent). 



Ground State, Probability Density 



Figure 16: Ground state density as a function of site number. The probability 
density is obtained by multiplying the values in the ordinate by the individual 
resonant wavefunctions, which are considered to be highly peaked at each site. 
The density has a gaussian envelope and does not exhibit nodes. 
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Figure 17: Vectors describing a 2D array. The components of b^ and are 
given in the text. 



generated by ai = (V3,0),a2 = (- ^3/2, 3/2), ag = (-^/3/2, -3/2) (grid A) 



while the other sublattice is obtained by adding the vectors bi = (0,1), b2 = 
(-y3/2,-l/2),b3 = (V3/2,-l/2). These vectors are so far dimensionless. 
The position vectors r^,rB of the periodic lattices are obtained by introducing 
the factor A (with the dimensions of length). Deformed lattices can also be 
described by these vectors, but the position vectors become more complicated 
functions of a^ , b^ . We denote by A the vector parametrizing sublattice A. For 
B we use A + bi. The state vectors (eigenvectors) for individual potential wells 
on grid A shall be denoted by |A), giving wave functions of individual wells as 
^yi(r — r^) = (r|A). For grid B we use |A + bi). The tight binding hamiltonian 
in this case is given by 



H 



a^|A)(A| +/3^|A + bi)(A + bi| + 



A A 



+ 




(121) 



A, 1=1,2,3 



The usual Pauli operators are constructed through the definitions 




CT- 



cr. 



.t 

+ 



(122) 



A 




(123) 



A 



while the kinetic operators 11, 11^ are defined as 
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VERTICAL SIDES VERTICAL SIDES 

ARE NOT DEFORMED ^RE NOT DEFORMED 




REGULAR HEXAGON REGULAR HEXAGON 

AS IWrriAL CELL AS INITIAL CELL 



Figure 18: Two possible deformations ol the lattice 



n = ^ A (I A) (A + b, - bi I + I A + bi) (A + b,|) . (124) 

A,i 

The spectrum and eigenfunctions are obtained again by squaring H. With M 
and Eq given as before, we obtain 



H = Eo + Ma3 + a+U + a^U'' (125) 

and 

{H - Eo f = + ntf (126) 
The spectrum and eigenfunctions are then 



E(k):^Eo± /A2|^e»2'^^b..k|2 + ^/2 (i27) 



It is well known that the degeneracy points of the spectrum for the massless 
case are kg = ±^(1, — a/S)- Around such points one finds 



E{k - ko) -Eo = ±v/A2fc2 + Af2 (^29) 
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^^^^^ 




Figure 19: Resonators with high dieletric constant e ~ 34. Courtesy of F. 
Mortessagne 

5.2.2 Tight binding and approximate isotropy 

We claim that rotational symmetry around conical points is a direct conse- 
quence of the tight binding approximation, as we shall see. It is well known 
that rotational symmetry in the Dirac equation demands a transformation of 
both orbital and spinorial degrees of freedom. It is in the orbital part that 
we shall concentrate by studying the energy surfaces around degeneracy points 
beyond the tight binding model. In our study, it will suffice to look inside the 
first Brillouin zone since the rest of the reciprocal lattice can be obtained by 
periodicity. Small deviations from degeneracy points (denoted by ko) in the 
form k = ko + K give the energy 

E = A|^exp(a(ko + k) -bi)! ~ AA|k;|, (130) 

i 

which is rotationally invariant in k. 

A second-neighbor interaction of strength A' modifies the kinetic operator 
n as 

n = A J2 ^b. +A' J2 ^a. +T_a., (131) 
i=l,2,3 1=1,2,3 

where the vectors a^ have now appeared, connecting a point with its six second 
neighbors. The energy equation becomes 

£; = I A ^ exp (iAk • b^) -I- A' ^ 2 cos (Ak • a,) | . (132) 

i i 

We expect a deviation of degeneracy points kp, for which k = kg -)- k. Upon 
linearization of the exponentials in k we find the energy 
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Figure 20: Energy surfaces computed form our dispersion formula in a periodic 
lattice. The conical points are visible at the six corners of the Brillouin zone. 



E ~ V(k • u)2 + (/« . v)2 (133) 
where the vectors are given by 

u = AA^cos(Ak;)-b,)b, (134) 

i 

V = AA ^ sin(Ak[, • bi)b, + 2AA' ^ sin(Ak;, • a,)a, (135) 

i i 

The presence of A' gives the energy surfaces (|133p as cones with elliptic 
sections whenever n is inside the first Brillouin zone. Regardless of how we 
complete the energy contours to recover periodicity, it is evident that the re- 
sulting surfaces are not invariant under rotations around degeneracy points. 
The circular case is recovered only when A' = 0, leading to kp = ko. In this 
case, the vectors reduce to v = (1, 0), u = (0, 1) when ko is the degeneracy point 
at (1/2A,0). 

In summary, extending the interactions to second neighbors has the effect 
of breaking the isotropy of space AROUND CONICAL POINTS, which is an 
essential property of the free Dirac theory. 

5.3 The Dirac oscillator in 2D 

We deform the lattice through an extension of the kinetic operators, just as 
in the one dimensional case. Let us consider site dependent transition ampli- 
tudes A(A, A + bi) connecting the sites labeled by A, A + bi. Again, these 
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Figure 21: First neighbor interaction, circular contours near the corners of the 
first Brillouin zone 




Figure 22: Second neighbor interaction, elhptic contours 
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are related to distances c?(A, A + bi) between resonators as A(A, A + bi) = 
Aexp(— (i(A, A + bi)/A). Now we define the ladder operator 

ar^Y. ^(A, A + b,) (|A)(A + b, - bil + |A + bi)(A + b,|) (136) 

A,i 

and impose [ar,aj] = u. After some algebra, one can prove that this leads to 
three recurrence relations. The first relation is 

A(A,A + bi)^A, (137) 

meaning that the vertical distances are fixed as a constant (the coupling is a 
constant A). The second and third relations give 

A2(A,A + b2) + A2(A + b2,A + b2-b3) = (138) 
A2(A + bi, a + bi + ba) + A2(A + bz + bi, A + bi + ba - bg), 



A2(A,A + b2) + A2(A,A + b3) = (139) 
A^(A + bi, A + bi - ba) + A^(A + bi, A + bi - ba) + 

It is the third relation what gives the scaling of distances in terms of our fre- 
quency uj: distances should increase in order to satisfy (|140p . The second re- 
lation simply establishes the equality of the lengths of opposite sides for each 
hexagonal cell. These relations seem to be complicated, but one can use a pro- 
gram to generate all lattice points consistently. We do so by starting with a 
regular hexagonal cell. The analogy between our model and the 2 dimensional 
Dirac oscillator becomes exact when the number of resonators is large. 




b 10 15 20 25 30 



Figure 23: A lattice produced with our recurrence relation. A regular hexagonal 
cell is used as a seed. A resonator is placed on each vertex of the lattice. A 
choice of deformation angle may produce periodicity in one direction (trivial 
case) 
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5 10 15 20 25 SO 



Figure 24: A lattice produced with our recurrence relation. Resonators are 
placed at the vertices of the array. A regular hexagonal cell is used as a seed 
(at the origin). No periodicity. 

The resulting hamiltonian of this problem is 

H = Eo + (T3M + (T+ar + (7-4 (140) 

with eigenvalues 



E{Nr) = Eo± ^/uj{Nr + l) + M^, 0<Nr< A'^/uj (141) 

The results obtained so far confirm our suspicion that the spectrum around 
conical points becomes more spaced with a square root law. See the figures for 
reflection and transmission measurements between antennas in the array. The 
peaks are localized around the blue cone located at the resonance, where the 
Dirac point should lie. 

Summarizing the results of this section, we have formulated a Dirac equa- 
tion in hexagonal lattices and justified the use of tight binding arrays, together 
with an experimental evidence of nearest neighbor coupling and its exponential 
decay as a function of separation distance. We provided a useful description 
for a problem motivated by graphene and the emulation of Dirac-Moshinsky 
oscillators in electromagnetic billiards. Moreover, we have developed a method 
to analyze deformations through the algebraic properties of the system, an idea 
that opens a window for the realization of other integrable systems. The ex- 
perimental realization of the DM0 depends crucially on the measured reflection 
peaks, as shown in the preliminary experimental results in the figure. So far, 
the location of Dirac points has been successful [55] and the distortion of the 
spectrum upon deformations is also visible. It is left to run more experiments 
in order to have a clear indication of a square root law for the spectrum and a 
localization of wavefunctions provided by the constantly increasing distance be- 
tween resonators. It must be mentioned that a large number of such resonators 
in our setup is mandatory, allowing the possibility of neglecting finite-size and 
boundary effects. 
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Graphen (easel) S12 




Figure 25: Preliminary experimental results for the transmission between an- 
tennas in the array of figure [221 as a function of the frequency (GHz). The 
blue line indicates the Dirac point. The equally spaced spectrum appears due 
to the deformation. The gap indicates the zero point energy of the oscillator. 
Courtesy of F. Mortessagne. 




Figure 26: Preliminary experimental results for the reflection in the array of 
figure [23l as a function of the frequency (GHz). The spaced spectrum appears 
due to the deformation (similar to a square root law, but with asymmetries). A 
blue curve indicates the location of the Dirac point. The gap indicates the zero 
point energy of the oscillator. Courtesy of F. Mortessagne 
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6 Conclusion 



We have reviewed the subject starting from its simple formulation as a potential 

problem, then revisited some of its achievements in hadron spectroscopy and fi- 
nally gave the reasons why this system can be realized in nature by careful 
constructions. The original system proposed by Moshinsky is simple enough to 
be considered as a paradigmatic model. Yet it possesses a richness of interpre- 
tations provided by its many formulations (in many dimensions and for many 
particles), by its original applicability to bound and composite systems and, in 
recent times, by the analogies that can be established in connection with two ar- 
eas that are active and prolific: two dimensional materials and quantum-optical 
traps. 
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